We wish to create conservation-sensitive, climate-smart fisheries closures for the ABNJ of the Pacific Ocean using prioritizr.
Prioritizr (Hanson et al., 2018) is a spatial prioritization software (R package) that uses integer linear programming (ILP) allowing it to create exact solutions to large spatial planning problems faster than its heuristic counterparts (e.g. Marxan).
This file runs multiple scripts found in the main Github repository: tinbuenafe/PacificProj, where the codes are explained in greater detail.
Before running prioritizr, we first had to assemble the following data layers:
The study area has the following characteristics:
We ran the following code using polygons from rnaturalearth package and the EEZ .shp files from Marine Regions (Flanders Marine Institute, 2019) to create a Pacific-centered map of the global ABNJs.
The relevant scripts for creating the study area are labeled with a prefix of 01_StudyArea in ~/scripts/.
# to create the Pacific-centered maps:
source("scripts/01a_StudyArea_CreatingPUs.R")
# this script uses the following functions (embedded actual script):
# function to create polygon of the ABNJ by masking the EEZs and the landmasses
source("scripts/study_area/fCreateMaskedPolygon.R")
# function to reproject these polygon to Robinson projection
source("scripts/study_area/fConvert2PacificRobinson.R")
# function to create the hexagonal planning units
source("scripts/study_area/fCreate_PlanningUnits.R")
This figure shows the Pacific-centered global map (using Robinson projection).
world_robinson <- readRDS("outputs/01_StudyArea/01a_StudyArea/PacificCenterLand.rds")
global_plot <- ggplot() +
geom_sf(data = world_robinson, colour = "grey20", fill="grey20", size = 0.1, show.legend = "line") +
theme_bw()
global_plot
We then created equal-area hexagonal planning units, or PUs, (0.5\(^\circ\) x 0.5\(^\circ\) at the equator) of the study area. This resulted in 31,917 PUs.
# define objects for plotting; objects are the result of running 01_StudyArea
pacific_robinson <- readRDS("outputs/01_StudyArea/01a_StudyArea/PacificCenterABNJ.rds")
PUsPac <- readRDS("outputs/01_StudyArea/01a_StudyArea/PacificABNJGrid_05deg.rds")
Bndry <- fCreateRobinsonBoundary(west = 78, east = 140, north = 51, south = 60) # coordinates in degrees
study_area <- ggplot() +
geom_sf(data = pacific_robinson, colour = NA, fill = NA, size = 0.2, show.legend = "line") +
geom_sf(data = world_robinson, color = "grey20", fill="grey20", size = 0.1, show.legend = "line") +
geom_sf(data = PUsPac, colour = "grey64", aes(fill = "ABNJ"), size = 0.1, show.legend = TRUE) +
scale_fill_manual(name = "Study Area", values = c("ABNJ" = "coral3")) +
coord_sf(xlim = c(st_bbox(Bndry)$xmin, st_bbox(Bndry)$xmax), # Set limits based on Bndry bbox
ylim = c(st_bbox(Bndry)$ymin, st_bbox(Bndry)$ymax),
expand = TRUE) +
theme_bw()
study_area
The following were used in this study:
First, we compiled a list of bycatch reported in the Pacific ABNJ. Then, we used the IUCN distribution maps for these species. These polygons are from the IUCN Red List, showing the known ranges of the species (IUCN Red List, 2021).
We limited the bycatch included in this study to sea turtles, but this workflow can be easily replicated and expanded to include more species.
The following were the sea turtles reported as bycatch and intersected with the study area:
The relevant scripts for creating the bycatch distribution layer are labeled with a prefix of 04_IUCN in ~/scripts/.
The plots below show the distribution of these species using fIUCNINtersect().
# the function is found in:
source("scripts/04a_IUCN_fIUCNINtersect.R")
# the function is ran in:
source("scripts/04b_IUCN_fIUCNIntersectRun.R")
Just like the bycatch layer, we compiled a list of commercially targeted large pelagic species managed in the Pacific from different sources. From this list, we have 4 species with available data on their global larval spawning grounds:
The models for these species are from a report by Mercer (2019) where he created Generalized Additive Models (GAMs) to predict spawning grounds using (Nishikawa et al., 1985). He used mean seasonal abundance from a historical dataset of plankton tows of fish larvae (dated back to the 1950s), and oceanographic and environmental conditions from global databases.
We reran his models (see scripts labeled with a prefix 05a_Commercial in ~/scripts/ for complete reruns) and used the medians of the predictions for each species as probability thresholds to create the distribution maps and the function fCommercialFeat().
# to intersect the results of the GAMs with the PUs
# the function is in:
source("scripts/05b_Commercial_fCommercialFeat.R")
# the function is ran at:
source("scripts/05c_Commercial_fCommercialFeatRun.R")
global_plots
This is an example of how the data looks like (so far) for bycatch:
head(bycatch)
## Simple feature collection with 6 features and 2 fields
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: -3779074 ymin: 168708.9 xmax: -3695824 ymax: 473118.3
## CRS: +proj=robin +lon_0=180 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs
## # A tibble: 6 x 3
## cellsID geometry species
## <int> <POLYGON [m]> <chr>
## 1 1 ((-3751324 216773.5, -3779074 232795.1, -3779074 26483… Caretta_caret…
## 2 2 ((-3751324 312902.8, -3779074 328924.4, -3779074 36096… Caretta_caret…
## 3 3 ((-3751324 409032.1, -3779074 425053.6, -3779074 45709… Caretta_caret…
## 4 4 ((-3723574 168708.9, -3751324 184730.5, -3751324 21677… Caretta_caret…
## 5 5 ((-3723574 264838.2, -3751324 280859.7, -3751324 31290… Caretta_caret…
## 6 6 ((-3723574 360967.5, -3751324 376989, -3751324 409032.… Caretta_caret…
And for the large pelagics:
head(pelagic)
## Simple feature collection with 6 features and 2 fields
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: -3779074 ymin: 168708.9 xmax: -3695824 ymax: 473118.3
## CRS: +proj=robin +lon_0=180 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs
## # A tibble: 6 x 3
## cellsID geometry species
## <int> <POLYGON [m]> <chr>
## 1 1 ((-3751324 216773.5, -3779074 232795.1, -3779074 264838.2, -3… ALB
## 2 2 ((-3751324 312902.8, -3779074 328924.4, -3779074 360967.5, -3… ALB
## 3 3 ((-3751324 409032.1, -3779074 425053.6, -3779074 457096.7, -3… ALB
## 4 4 ((-3723574 168708.9, -3751324 184730.5, -3751324 216773.5, -3… ALB
## 5 5 ((-3723574 264838.2, -3751324 280859.7, -3751324 312902.8, -3… ALB
## 6 6 ((-3723574 360967.5, -3751324 376989, -3751324 409032.1, -372… ALB
Each sf object has the following columns:
cellsID is the unique ID for each hexagonal PUspecies is the species code (e.g. ALB = Albacore)geometryTo create climate-smart closures, we prioritized protection of areas with:
We determined these areas using the following metrics:
Relative Climate Exposure (RCE) index
\(\large RCE(yr^{-1}) = \frac{{\text Slope \hspace{2mm}} (^{\circ}Cyr^{-1} \hspace{2mm} 2050 - 2100)}{{\text Seasonal \hspace{2mm} Range} \hspace{2mm} (^{\circ}C \hspace{2mm} 2015 - 2020)}\)
Areas with low RCE index values are areas that will be exposed to low levels of warming relative to their current seasonal fluctuations of temperature. Areas of low RCE index values are prioritized for selection.
Climate velocity
\(\large Climate \hspace{2mm} velocity \hspace{2mm} (kmyr^{-1}) = \frac{{\text Slope \hspace{2mm}} (^{\circ}Cyr^{-1} \hspace{2mm} 2050 - 2100)}{{\text Spatial \hspace{2mm} gradient} \hspace{2mm} (^{\circ}Ckm^{-1} \hspace{2mm} 2050 - 2100)}\)
Areas of slow climate velocity are areas that are more likely to retain their current environmental conditions and consequentially, their biodiversity. Areas of slow climate velocity are prioritized for selection.
The temperatures used to calculate the aforementioned metrics were derived from 12 Global Circulation Models (GCMs) from the Couple Model Intercomparison Project Phase 6 (CMIP6; https://esgf-node.llnl.gov). The GCMs are under these three climate scenarios:
We intersected the RCE values with the PUs using fClimateInt(). The relevant scripts for creating the climate-smart features layer are found in ~/scripts/ with the prefix 07_Climate.
# function is found at:
source("scripts/07a_Climate_fClimateInt.R")
# function is ran at:
source("scripts/07b_Climate_fClimateIntRun.R")
We intersected the climate velocity values with the PUs using the same function, fClimateInt().
The cost layer represents the opportunity cost to fishing of establishing the proposed fisheries closures.
The cost layer for each planning unit was represented in US Dollars, and was calculated by multiplying the:
In addition, we smoothed the cost layer using focal() to reduce the splotchiness of the data.
# function is found at:
source("scripts/06a_Cost_fCostPU.R")
# function is ran at:
source("scripts/06b_Cost_fCostPURun.R")
## [1] "minimum: 0.00668002245947719 maximum: 813.370788574219"
We used log-transformed cost for the figure below:
We join the climate-smart layer with the conservation features layer. Then, we only retain planning units that have RCE and/or climate-velocity values belonging to their lower 25th percentile range. This was how we priortized protection of areas that are predicted to have lower exposure to climate warming and higher retention of biodiversity.
We used the function fFilterQuartile() and the scripts that are relevant for this step are found in ~/scripts/ with the prefix 08_Filter.
# function is found at:
source("scripts/08a_Filter_fFilterQuartile.R")
# function is ran at:
source("scripts/08d_Filter_fFilterQuartile25Run_NoProv.R")
This was done for both the bycatch and large pelagic layers, across different scenarios (SSP1-2.6, SSP2-4.5, SSP5-8.5). Below is an example of how the data looks like, with the following columns:
cellsID: unique ID for each PUspecies: species codetotal_area: the total area covered by each subtyped feature (in km\(^{-1}\))velocity: climate velocityvelo_tvalue: transformed values of climate velocity (velocity * 10)RCE: RCE index valuesRCE_tvalue: transformed values of RCE (cube root of RCE)# for bycatch under the climate scenario SSP1-2.6.
bycatchSSP126 <- readRDS("outputs/08_Filter/08d_Filter25_NoProv/bycatchIUCNSSP126_25percentile.rds") %>%
as_tibble() %>%
dplyr::select(-area_km2, -cellsID.2, -geometry)
head(bycatchSSP126)
## # A tibble: 6 x 7
## cellsID new_features velocity velo_tvalue RCE RCE_tvalue total_area
## <int> <chr> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 12 Caretta_caretta_IU… 0.239 2.39 0.0352 0.328 19761581.
## 2 21 Caretta_caretta_IU… 0.143 1.43 0.0239 0.288 19761581.
## 3 22 Caretta_caretta_IU… 0.0163 0.163 0.00493 0.170 19761581.
## 4 23 Caretta_caretta_IU… -0.0532 -0.532 0.0107 0.221 19761581.
## 5 32 Caretta_caretta_IU… 0.237 2.37 0.0355 0.329 19761581.
## 6 33 Caretta_caretta_IU… 0.0717 0.717 0.0131 0.236 19761581.
# for large pelagics under the cliamte scenario SSP5-8.5.
commercialSSP585 <- readRDS("outputs/08_Filter/08d_Filter25_NoProv/commercialSSP585_25percentile.rds") %>%
as_tibble() %>%
dplyr::select(-area_km2, -cellsID.2, -geometry)
head(commercialSSP585)
## # A tibble: 6 x 7
## cellsID new_features velocity velo_tvalue RCE RCE_tvalue total_area
## <int> <chr> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 84 ALB 37.1 371. 2.24 1.31 5127127.
## 2 86 ALB 38.7 387. 1.18 1.06 5127127.
## 3 88 ALB 39.1 391. 2.19 1.30 5127127.
## 4 89 ALB 46.5 465. 1.29 1.09 5127127.
## 5 94 ALB 62.2 622. 1.40 1.12 5127127.
## 6 95 ALB 42.7 427. 1.19 1.06 5127127.
Aside from these climate-smart closures, we also explored climate-uninformed, wherein the climate-smart features were not included in the analyses. Subsequently, the step of filtering and retention of planning units was also not included.
source("scripts/08e_Filter_fFilterQuartile100Run_NoProv.R")
This is how the data looks like:
# for bycatch:
bycatch_uninformed <- readRDS("outputs/08_Filter/08e_Filter100_NoProv/bycatchIUCNuninformed_100percentile.rds") %>%
as_tibble() %>%
dplyr::select(-area_km2, -geometry)
head(bycatch_uninformed)
## # A tibble: 6 x 3
## cellsID new_features total_area
## <int> <chr> <dbl>
## 1 1 Caretta_caretta_IUCN 62747287.
## 2 2 Caretta_caretta_IUCN 62747287.
## 3 3 Caretta_caretta_IUCN 62747287.
## 4 4 Caretta_caretta_IUCN 62747287.
## 5 5 Caretta_caretta_IUCN 62747287.
## 6 6 Caretta_caretta_IUCN 62747287.
# for large pelagics:
commercial_uninformed <- readRDS("outputs/08_Filter/08e_Filter100_NoProv/commercialuninformed_100percentile.rds") %>%
as_tibble() %>%
dplyr::select(-area_km2, -geometry)
head(commercial_uninformed)
## # A tibble: 6 x 3
## cellsID new_features total_area
## <int> <chr> <dbl>
## 1 1 ALB 11708096.
## 2 2 ALB 11708096.
## 3 3 ALB 11708096.
## 4 4 ALB 11708096.
## 5 5 ALB 11708096.
## 6 6 ALB 11708096.
To determine the relative targets per subtyped feature, we checked each one’s IUCN category (IUCN, 2015) using rredlist(). Threatened features (i.e. classified as vulnerable, endangered, threatened) were assigned a fixed maximum target. Features (\(i\)) that were not classified as threatened (e.g. classified as near-threatened, etc.) were assigned a relative target based on their distribution.
Relative target\(_i\) \(= Target_{max} - \frac{PU_i}{PU_{total}} \cdot (Target_{max} - Target_{min})\)
PU_i: number of planning units of feature \(i\)PU_{total}: total number of planning units in planning regionTarget_{max}: Maximum target / the fixed maximum target assigned to threatened features (e.g. 1)Target_{min}: Minimum target (e.g. 0.1)Lower targets are assigned to features with broader distribution, and higher targets to features with more restricted distribution.
We reran this step for an array of targets, for both climate-smart and climate-uninformed closures. We multiplied 0.25 to each of the maximum representative targets (%) to determine the percentage of area, also known as the effective target (%) that is actually conserved after taking only retaining PUs less than the 25th percentile of both climate metrics.
| Runs | Maximum representative target (%) | Effective target (%) |
|---|---|---|
| 1 | 100 | 25 |
| 2 | 90 | 22.5 |
| 3 | 80 | 20 |
| 4 | 70 | 17.5 |
| 5 | 60 | 15 |
| 6 | 50 | 12.5 |
| 7 | 40 | 10 |
| 8 | 30 | 7.5 |
| 9 | 20 | 5 |
| 10 | 10 | 2.5 |
We used the function fRepresentTarget(), and the scripts are found in ~/scripts/ with the prefix 09.
# the function is found in:
source("scripts/09a_Target_fRepresentTarget.R")
# the runs for climate-smart closures are found in:
source("scripts/09d_Target_fRepresentTargetSmartRun_NoProv.R")
# the runs for climate-uninformed closures are found in:
source("scripts/09e_Target_fRepresentTargetUninformedRun_NoProv.R")
For this .rmd, we’ll be showing results using Runs 2 (Max = 0.9, Min = 0.1), 5 (Max = 0.6, Min = 0.1) and 8 (Max = 0.3, Min = 0.1). This workflow shows that it can be replicated using other parameters.
Below is an example of how the data looks like:
target90_commercial <- readRDS("outputs/09_Target/09d-e_TargetRuns/02_Target90/SSP126/target_commercialSSP126.rds") %>%
as_tibble()
head(target90_commercial)
## # A tibble: 4 x 5
## new_features total_area category geometry target
## <chr> <dbl> <chr> <MULTIPOLYGON [m]> <dbl>
## 1 ALB 3964054. NT (((-2058566 -3772591, -2086316 -37565… 86.3
## 2 SKP 2392837. LC (((-1559064 -2715169, -1586814 -26991… 87.8
## 3 SWO 4428216 LC (((-2058566 -3772591, -2086316 -37565… 85.8
## 4 YFT 3905366. NT (((-671059.4 -696454.6, -698809.5 -68… 86.3
target90_bycatch <- readRDS("outputs/09_Target/09d-e_TargetRuns/02_Target90/SSP126/target_bycatchSSP126.rds") %>%
as_tibble()
head(target90_bycatch)
## # A tibble: 5 x 5
## new_features total_area category geometry target
## <chr> <dbl> <chr> <MULTIPOLYGON [m]> <dbl>
## 1 Rep-1347 56020. VU (((-1753314 -2859363, -1781065 -28433… 90
## 2 Rep-2226 53352 EN (((-1753314 -2859363, -1781065 -28433… 90
## 3 Rep-3437 20644556. VU (((688697.1 -5743241, 660946.9 -57272… 90
## 4 Rep-4127 58687. CR (((-1753314 -2859363, -1781065 -28433… 90
## 5 Rep-5414 16006. VU (((-1864315 -2763234, -1892065 -27472… 90
new_features: featurestotal_area: the total area covered by each subtyped feature (in km\(^{-1}\))category: IUCN redlist categorytarget: assigned target in percentageprioritizr (v. 7.0.1.2) uses the following functions to define the problem.
# installing latest version of prioritizr from GitHub
devtools::install_github("prioritizr/prioritizr")
# installing gurobi as the optimizer
install.packages("/Library/gurobi911/mac64/R/gurobi_9.1-1_R_4.0.2.tgz", repos = NULL)
library(priortizr)
library(gurobi)
p1 <- prioritizr::problem(x, feature_list, "cost") %>%
add_min_set_objective() %>%
add_relative_targets(target_list) %>%
add_binary_decisions() %>%
add_gurobi_solver(gap = 0.1, verbose = FALSE)
s1 <- solve(p1)
problem(): we define the datax: data frame having all the data layers (all planning units with their specific ID (cellsID), features data, cost)feature_list: list of the features in xcost: name of the column in x where the cost is foundadd_min_set_objective(): function that allows us to employ the minimum set objective function. This allows prioritizr to create spatial plans that will minimize the cost.add_relative_targets(): target_list list of assigned targets for each featureadd_binary_decisions(): prioritizr will only either select or not select a PUadd_gurobi_solver(): Gurobi was used as the optimizer to solve the problem (at a default 10%)solve(): prioritizr solves the problem p1We used the fPrioritizrRun() function to streamline the runs. fPrioritizrRun() arranges all the data and makes sure that they’re in the right format, then runs problem() and solve(). This function also saves the summaries each plan (cost, % area protected and number of representative features that achieved 30% or more protection).
The relevant scripts for this step is found in ~/scripts/ with the prefix 10_Prioritizr.
# the function is in:
source("scrpits/10a_Prioritizr_fPrioritizr.R")
# to run the function with climate-smart data
source("scripts/10h_Prioritizr_fPrioritizrIUCNSmart_NoProv.R")
# to run the function with climate-uninformed data
source("scripts/10i_Prioritizr_fPrioritizrIUCNUninformed_NoProv.R")
# to create no-regret closures with the climate-smart plans
To create no-regret plans, we used thre function fCreateNoRegret() which intersected the selected PUs in all three climate scenarios. The relevant scripts for this step is found in ~/scripts with the prefix 11_NoRegret.
# the function is in:
source("scripts/11_NoRegret_fCreateNoRegret.R")
# the runs are found in:
source("scripts/11e_NoRegret_fCreateNoRegretIUCNRuns_NoProv.R")
As mentioned in the previous section, we will only be showing results of Runs 2 (Max = 0.9, Min = 0.1), 5 (Max = 0.6, Min = 0.1) and 8 (Max = 0.3, Min = 0.1).
SSP1-2.6 run
SSP126_target90 <- readRDS("outputs/10_Prioritizr/10h-i_IUCNRuns_NoProv/02_Target90/solution_SSP126.rds")
cost is in USDperc_area_protected is the percentage of the Pacific ABNJ protectedprotected_PU is the number of features that have >= 30% (by area) protection## # A tibble: 0 x 12
## # … with 12 variables: target_scenarios <chr>, sum_area <dbl>,
## # percent_area <dbl>, num_pu <dbl>, total_cost <dbl>, median_velocity <dbl>,
## # median_tvelocity <dbl>, median_RCE <dbl>, median_RCEtvalue <dbl>,
## # represented_features <dbl>, target <dbl>, scenario <chr>
SSP2-4.5 run
SSP245_target90 <- readRDS("outputs/10_Prioritizr/10h-i_IUCNRuns_NoProv/02_Target90/solution_SSP245.rds")
cost is in USDperc_area_protected is the percentage of the Pacific ABNJ protectedprotected_PU is the number of features that have >= 30% (by area) protection## # A tibble: 0 x 12
## # … with 12 variables: target_scenarios <chr>, sum_area <dbl>,
## # percent_area <dbl>, num_pu <dbl>, total_cost <dbl>, median_velocity <dbl>,
## # median_tvelocity <dbl>, median_RCE <dbl>, median_RCEtvalue <dbl>,
## # represented_features <dbl>, target <dbl>, scenario <chr>
SSP5-8.5 run
SSP585_target90 <- readRDS("outputs/10_Prioritizr/10h-i_IUCNRuns_NoProv/02_Target90/solution_SSP585.rds")
cost is in USDperc_area_protected is the percentage of the Pacific ABNJ protectedprotected_PU is the number of features that have >= 30% (by area) protection## # A tibble: 0 x 12
## # … with 12 variables: target_scenarios <chr>, sum_area <dbl>,
## # percent_area <dbl>, num_pu <dbl>, total_cost <dbl>, median_velocity <dbl>,
## # median_tvelocity <dbl>, median_RCE <dbl>, median_RCEtvalue <dbl>,
## # represented_features <dbl>, target <dbl>, scenario <chr>
noregret_target90 <- readRDS("outputs/11_NoRegret/11e_IUCNRuns_NoProv/noregretclosures_Target90.rds")
uninformed_target90 <- readRDS("outputs/10_Prioritizr/10h-i_IUCNRuns_NoProv/02_Target90/solution_uninformed.rds")
cost is in USDperc_area_protected is the percentage of the Pacific ABNJ protectedprotected_PU is the number of features that have >= 30% (by area) protection## # A tibble: 1 x 12
## target_scenarios sum_area percent_area num_pu total_cost median_velocity
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Target90_summaryuninf… 1.68e7 19.8 6308 120539. 16827221.
## # … with 6 more variables: median_tvelocity <dbl>, median_RCE <dbl>,
## # median_RCEtvalue <dbl>, represented_features <dbl>, target <dbl>,
## # scenario <chr>
SSP1-2.6 run
SSP126_target60 <- readRDS("outputs/10_Prioritizr/10h-i_IUCNRuns_NoProv/05_Target60/solution_SSP126.rds")
cost is in USDperc_area_protected is the percentage of the Pacific ABNJ protectedprotected_PU is the number of features that have >= 30% (by area) protection## # A tibble: 0 x 12
## # … with 12 variables: target_scenarios <chr>, sum_area <dbl>,
## # percent_area <dbl>, num_pu <dbl>, total_cost <dbl>, median_velocity <dbl>,
## # median_tvelocity <dbl>, median_RCE <dbl>, median_RCEtvalue <dbl>,
## # represented_features <dbl>, target <dbl>, scenario <chr>
SSP2-4.5 run
SSP245_target60 <- readRDS("outputs/10_Prioritizr/10h-i_IUCNRuns_NoProv/05_Target60/solution_SSP245.rds")
cost is in USDperc_area_protected is the percentage of the Pacific ABNJ protectedprotected_PU is the number of features that have >= 30% (by area) protection## # A tibble: 0 x 12
## # … with 12 variables: target_scenarios <chr>, sum_area <dbl>,
## # percent_area <dbl>, num_pu <dbl>, total_cost <dbl>, median_velocity <dbl>,
## # median_tvelocity <dbl>, median_RCE <dbl>, median_RCEtvalue <dbl>,
## # represented_features <dbl>, target <dbl>, scenario <chr>
SSP5-8.5 run
SSP585_target60 <- readRDS("outputs/10_Prioritizr/10h-i_IUCNRuns_NoProv/05_Target60/solution_SSP585.rds")
cost is in USDperc_area_protected is the percentage of the Pacific ABNJ protectedprotected_PU is the number of features that have >= 30% (by area) protection## # A tibble: 0 x 12
## # … with 12 variables: target_scenarios <chr>, sum_area <dbl>,
## # percent_area <dbl>, num_pu <dbl>, total_cost <dbl>, median_velocity <dbl>,
## # median_tvelocity <dbl>, median_RCE <dbl>, median_RCEtvalue <dbl>,
## # represented_features <dbl>, target <dbl>, scenario <chr>
noregret_target60 <- readRDS("outputs/11_NoRegret/11e_IUCNRuns_NoProv/noregretclosures_Target60.rds")
uninformed_target60 <- readRDS("outputs/10_Prioritizr/10h-i_IUCNRuns_NoProv/05_Target60/solution_uninformed.rds")
cost is in USDperc_area_protected is the percentage of the Pacific ABNJ protectedprotected_PU is the number of features that have >= 30% (by area) protection## # A tibble: 1 x 12
## target_scenarios sum_area percent_area num_pu total_cost median_velocity
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Target60_summaryuninf… 11203920 13.2 4200 73545. 11203920
## # … with 6 more variables: median_tvelocity <dbl>, median_RCE <dbl>,
## # median_RCEtvalue <dbl>, represented_features <dbl>, target <dbl>,
## # scenario <chr>
SSP1-2.6 run
SSP126_target30 <- readRDS("outputs/10_Prioritizr/10h-i_IUCNRuns_NoProv/08_Target30/solution_SSP126.rds")
cost is in USDperc_area_protected is the percentage of the Pacific ABNJ protectedprotected_PU is the number of features that have >= 30% (by area) protection## # A tibble: 0 x 12
## # … with 12 variables: target_scenarios <chr>, sum_area <dbl>,
## # percent_area <dbl>, num_pu <dbl>, total_cost <dbl>, median_velocity <dbl>,
## # median_tvelocity <dbl>, median_RCE <dbl>, median_RCEtvalue <dbl>,
## # represented_features <dbl>, target <dbl>, scenario <chr>
SSP2-4.5 run
SSP245_target30 <- readRDS("outputs/10_Prioritizr/10h-i_IUCNRuns_NoProv/08_Target30/solution_SSP245.rds")
cost is in USDperc_area_protected is the percentage of the Pacific ABNJ protectedprotected_PU is the number of features that have >= 30% (by area) protection## # A tibble: 0 x 12
## # … with 12 variables: target_scenarios <chr>, sum_area <dbl>,
## # percent_area <dbl>, num_pu <dbl>, total_cost <dbl>, median_velocity <dbl>,
## # median_tvelocity <dbl>, median_RCE <dbl>, median_RCEtvalue <dbl>,
## # represented_features <dbl>, target <dbl>, scenario <chr>
SSP5-8.5 run
SSP585_target30 <- readRDS("outputs/10_Prioritizr/10h-i_IUCNRuns_NoProv/08_Target30/solution_SSP585.rds")
cost is in USDperc_area_protected is the percentage of the Pacific ABNJ protectedprotected_PU is the number of features that have >= 30% (by area) protection## # A tibble: 0 x 12
## # … with 12 variables: target_scenarios <chr>, sum_area <dbl>,
## # percent_area <dbl>, num_pu <dbl>, total_cost <dbl>, median_velocity <dbl>,
## # median_tvelocity <dbl>, median_RCE <dbl>, median_RCEtvalue <dbl>,
## # represented_features <dbl>, target <dbl>, scenario <chr>
noregret_target60 <- readRDS("outputs/11_NoRegret/11e_IUCNRuns_NoProv/noregretclosures_Target60.rds")
uninformed_target30 <- readRDS("outputs/10_Prioritizr/10h-i_IUCNRuns_NoProv/08_Target30/solution_uninformed.rds")
cost is in USDperc_area_protected is the percentage of the Pacific ABNJ protectedprotected_PU is the number of features that have >= 30% (by area) protection## # A tibble: 1 x 12
## target_scenarios sum_area percent_area num_pu total_cost median_velocity
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Target30_summaryuninf… 5601960 6.58 2100 32055. 5601960
## # … with 6 more variables: median_tvelocity <dbl>, median_RCE <dbl>,
## # median_RCEtvalue <dbl>, represented_features <dbl>, target <dbl>,
## # scenario <chr>
Since we ran prioritizr from 10 - 100% maximum target representation (in increments of 10%), we’re showing a summary of how cost and percent protection varied with increasing target representation.
source("scripts/12d_Stat_fCreateSummaryRuns.R")
# showing the log(cost), log(% area protected), number of protected planning units, and number of adequately represented features.
# this was done to determine which maximum representation targets to report/use in the runs
smart_plans
We also determined Cohen’s Kappa Coefficient between plans created for each run. Below are the summaries for each run as well as a matrix plot of the correlations between plans. We also show the summaries of the representation of each feature for each spatial plan.
The relevant scripts for this part are found in ~/scripts/ with prefix 12_Stat.
# the function for creating the correlation plots is found in:
source("scripts/12e_Stat_fKappaCorrplot.R")
# the runs are found in:
source("scripts/12f_Stat_fKappaCorrplotRun.R")
# the function for creating circular bar plots to show summaries of the representation of each feature:
source("scripts/12g_Stat_fCreateCircBarPlot.R")
# the runs are found in:
source("scripts/12h_Stat_fCreateCircBarPlotRuns.R")
Summary for Run 2: Maximum Target = 90%
## # A tibble: 4 x 12
## target_scenarios sum_area percent_area num_pu total_cost median_velocity
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Target90_summarySSP126 2.45e7 28.8 9203 324294. 0.0821
## 2 Target90_summarySSP245 2.98e7 35.0 11161 349721. 2.30
## 3 Target90_summarySSP585 3.45e7 40.5 12927 387218. 13.1
## 4 Target90_summaryuninf… 1.68e7 19.8 6308 120539. 16827221.
## # … with 6 more variables: median_tvelocity <dbl>, median_RCE <dbl>,
## # median_RCEtvalue <dbl>, represented_features <dbl>, target <dbl>,
## # scenario <chr>
Summary for Run 5: Maximum Target = 60%
## # A tibble: 4 x 12
## target_scenarios sum_area percent_area num_pu total_cost median_velocity
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Target60_summarySSP126 1.60e7 18.8 5994 168920. 0.0791
## 2 Target60_summarySSP245 1.81e7 21.2 6773 177037. 2.28
## 3 Target60_summarySSP585 2.09e7 24.5 7829 186072. 12.9
## 4 Target60_summaryuninf… 1.12e7 13.2 4200 73545. 11203920
## # … with 6 more variables: median_tvelocity <dbl>, median_RCE <dbl>,
## # median_RCEtvalue <dbl>, represented_features <dbl>, target <dbl>,
## # scenario <chr>
Summary for Run 8: Maximum Target = 30%
## # A tibble: 4 x 12
## target_scenarios sum_area percent_area num_pu total_cost median_velocity
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Target30_summarySSP126 7.74e6 9.09 2900 65808. 0.0768
## 2 Target30_summarySSP245 9.20e6 10.8 3450 71797. 2.36
## 3 Target30_summarySSP585 1.03e7 12.1 3872 73878. 12.9
## 4 Target30_summaryuninf… 5.60e6 6.58 2100 32055. 5601960
## # … with 6 more variables: median_tvelocity <dbl>, median_RCE <dbl>,
## # median_RCEtvalue <dbl>, represented_features <dbl>, target <dbl>,
## # scenario <chr>